#0  rarray.append(Rapameter("tau",  "Transcription-translation delay",    10.0,True,True,0.1));
#1  rarray.append(Rapameter("rp",   "Production rate (P)",                1.9,True,True,0.1));
#2  rarray.append(Rapameter("rm",   "Production rate (M)",                0.5,True,True,0.1));
#3  rarray.append(Rapameter("qp",   "Degradation rate (P)",               0.4,True,True,0.1));
#4  rarray.append(Rapameter("qm",   "Degradation rate (M)",               0.4,True,True,0.1));
#5  rarray.append(Rapameter("m",    "Non-linearity",                      3.5,True,True,0.1));
#6  rarray.append(Rapameter("n",    "Hill coefficient",                   6.0,True,True,2.1));
#7  rarray.append(Rapameter("k",    "Scaling constant",                   1.0,True,True,0.1));
from pylab        		import *
#from rapameters   		import *
from oscillator   		import *
from oscillator_hash  import *
from copy         		import *
from time         		import *
from raptools4    		import *
from wave         		import Wave
from random       		import *;

#TODO: ce ugotovis, da ne konvergira, spremeni resolucijo.
#TODO: na koncu so vsi toliko blizu, da se boljsi iscejo pocasneje. spremeniti korak?


if len(sys.argv) < 2:
  print "usage: rap.py amp wl"
  quit()

#CILJNI PARAMETRI
AMP = float(sys.argv[1])
WL  = float(sys.argv[2])
FRQ = 1.0/AMP;

#KORAK IN NAPAKA
STEP    = 0
AMP_MAX  = AMP*10
EPS_AMP = 0.05
EPS_FRQ = 0.01
EPS_TIM = 5

#OSCILATORJI
ROUNDING            = det_rounding(AMP);
SAMPLING_EXP        = 0.1
TIME                = int(math.ceil(WL*10));
DT                  = 1*10**(-SAMPLING_EXP);

#POPULACIJA
POPULATION          = [];
POP_MAX             = 400;
POP_LIFETIME        = 5;
POP_SIZE            = 50;
POP_SIZES           = []

#VERJETNOSTI
#               [tau, rp,   rm, qp, qm, m,  n,  k]
RECOMBINATION   = 1.0/8.0;
MUTATION_MASK   = [1,1,1,1,1,0,0,0]
MUTATION        = 1.0/5.0;


#ZACETNE VREDNOSTI
#               [tau, rp,   rm, qp, qm, m,  n,  k]
#rarray        = [float(WL)/2.0,1.0,1.0,1.0,1.0,3.0,2.0,1.0];
rarray        = [4.0,1.0,1.0,0.21,0.21,3.0,2.0,1.0];
rarray        = [1.0,1.0,1.0,1.0,1.0,3.0,2.0,1.0];


#ZAGON OSCILATORJA
start         = Oscillator(rarray,TIME,DT,ROUNDING, AMP_MAX)
WAVE          = Wave(float(AMP),1.0/float(WL),TIME,DT)
POPULATION    = [ Oscillator( gen_osc_repr(start, MUTATION, MUTATION_MASK), TIME, DT, ROUNDING, AMP_MAX) for i in range (POP_SIZE)];

for oscp in POPULATION:
  print oscp

##################
#ZACETEK RACUNANJA
##################
print "Time: ", TIME, " h"
print "Resolution: ", DT, "h"

while True:

  DT                  = 1*10**(-SAMPLING_EXP);


  print "-------------------------------------------------------------"
  print "STEP: ", STEP, " Time: ", TIME, " h", " Resolution: ", DT, "h"
  print "-------------------------------------------------------------"

  if STEP == 0:
    MUTATION_MASK   = [1,1,1,1,1,1,1,1]
    FITNESS_MASK    = [1,0,0]
    MUTATION        = 1.0/8.0
  elif STEP == 1:
    MUTATION_MASK   = [1,1,1,1,1,1,1,1]
    FITNESS_MASK    = [0,1,0]
    MUTATION        = 1.0/8.0
    DT  = 1*10**(-SAMPLING_EXP);
  elif STEP == 2:
    MUTATION_MASK   = [1,1,1,1,1,1,1,1]
    FITNESS_MASK    = [1,1,1]
    MUTATION        = 1.0/8.0  
  else:
      break;

  MUTATION        = 1.0/sum(MUTATION_MASK);

  #zacetni oscilator
  #figure(1)
  #spt = start.get_times();
  #spv = start.get_vals();
  #plot(spt, spv, linewidth=1.0, color = 'blue')
  #show();
  #quit()

  cmnbest     = float("inf");
  cmnworst    = 0;
  oscbest     = None
  oscworst    = None
  cut         = 0;


  #Zanka skozi generacije
  for GEN in range(POP_LIFETIME):

    print "#STEP: ",STEP," GENERATION (",str(GEN),")#, size: ", len(POPULATION)
    NEW_POPULATION = []
    POP_SIZES = POP_SIZES + [len(POPULATION)]

    
    #######################
    #KRITERIJ ZA PREZIVETJE
    #######################

    mn_list = []
    for osc in POPULATION:
      c = cmp_wav_osc_fitness(WAVE, osc, FITNESS_MASK)
      if c < float("inf"):
        mn_list = mn_list + [c];
      cut = mean(mn_list);  
    
    #Prepolovitev populacije
    if len(POPULATION) > POP_MAX:
      cut = median(mn_list)
      
    ##########
    #SELEKCIJA
    ##########
    cmnbest     = float("inf");
    cmnworst    = 0;
    for osc in POPULATION:

      cmn = cmp_wav_osc_fitness(WAVE, osc, FITNESS_MASK);

      #print id(osc), str(osc)
      #print id(osc), 'compare ', cmp_wav_osc_lst(WAVE, osc), 'mean', cmn
      #print id(osc), 'makes the cut(',str(cut),')? ', cmn <= cut;
      #print ''

      ####
      #CUT
      ####
      '''
      if cmn <= cut:
        if cmn < cmnbest:
          cmnbest = cmn;
          oscbest = osc;
        if cmn >= cmnworst:
          cmnworst = cmn;
        
        NEW_POPULATION = NEW_POPULATION + [osc]
        if size(NEW_POPULATION) > POP_MAX:
          print 'remove'
          NEW_POPULATION.remove(cmnworst)
      '''
      
      '''
      if cmn < cmnbest:
          cmnbest = cmn;
          oscbest = osc;
      if cmn >= cmnworst:
          cmnworst = cmn;
          oscworst = osc;
        
      NEW_POPULATION = NEW_POPULATION + [osc]
      if size(NEW_POPULATION) > POP_MAX:
          NEW_POPULATION.remove(oscworst)
          for onp in NEW_POPULATION:
              cmnworst = 0
              cmnnp = cmp_wav_osc_fitness(WAVE, onp);
              if cmn >= cmnworst:
                  cmnworst = cmnnp;
                  oscworst = onp;
      '''

      #Ce oscilator prezivi rez, ga dodaj v populacijo.
      #Doloci najboljsi in najslabsi primerek.
      if cmn <= cut:
          if cmn < cmnbest:
              cmnbest = cmn;
              oscbest = osc;
          if cmn >= cmnworst:
              cmnworst = cmn;
              oscworst = osc;
        
          #Ce velikost populacije presega maksimum, odstrani najslabsi primerek. 
          #Doloci nov najslabsi primerek.
          #TODO: Bottleneck.
          NEW_POPULATION = NEW_POPULATION + [osc]
          if size(NEW_POPULATION) > POP_MAX:
              NEW_POPULATION.remove(oscworst)
              for onp in NEW_POPULATION:
                  cmnworst = 0
                  cmnnp = cmp_wav_osc_fitness(WAVE, onp, FITNESS_MASK);
                  if cmn >= cmnworst:
                      cmnworst = cmnnp;
                      oscworst = onp;
      #########
      #DBGPLOT#
      #########
      #t = oscbest.get_times();
      #v = oscbest.get_vals();
      #plot(t,v,color=str(clr))
    
    #print 'before repr: ', len(NEW_POPULATION)

    ##############################
    #POPRAVLJANJE MUTACIJSKE MASKE
    ##############################  
    '''  
    if math.fabs(oscbest.get_amp()-float(AMP))/float(AMP)  <= MUTATION_CUT:
      print 'MUTATION CUT'
      MUTATION_MASK = [1,0,0,0,0,0,0,0]
      MUTATION = 1.0
    else:
      MUTATION_MASK = [1,1,1,1,1,1,1,1]
      MUTATION = 1.0/8.0
    '''
    
    #########################
    #REPRODUKCIJA IN MUTACIJE
    #########################
    off_size = 1+len(NEW_POPULATION);
    fitness_isum = sum([1.0/cmp_wav_osc_fitness(WAVE, osc, FITNESS_MASK) for osc in NEW_POPULATION])
    off_list = [[id(osc), osc, math.floor((1.0/cmp_wav_osc_fitness(WAVE, osc, FITNESS_MASK))/fitness_isum*off_size)] for osc in NEW_POPULATION]
    
    #print off_list
    #quit()
    
    for osc_repr in off_list:
      try:
        for i in range(int(osc_repr[2])):
          NEW_POPULATION = NEW_POPULATION + [ Oscillator( gen_osc_repr(osc_repr[1], MUTATION, MUTATION_MASK), TIME, DT, ROUNDING, AMP_MAX) ]
      except ValueError:
        print osc_repr[1]
        print osc_repr[2]
        continue
    
    #print off_list
    #print
    #print
    #print
    #print NEW_POPULATION
    POPULATION = NEW_POPULATION;

    ###########################
    #KONEC RACUNANJA GENERACIJE
    ###########################
  
  #################
  #NAJBOLJSI OSEBEK
  #################
  print "----------------------------------------------------------"
  print 'WAVE (RED)'
  print "----------------------------------------------------------"
  print 'amp: ',WAVE.get_amp();
  print 'frq: ', WAVE.get_frq();
  print 'time up: ',WAVE.get_tim_up();
  print 'time down: ',WAVE.get_tim_down();
  print "----------------------------------------------------------"
  print 'OSCILLATOR : CURRENT BEST', id(oscbest), ' (BLUE)'
  print "----------------------------------------------------------"
  print 'amp: ',oscbest.get_amp();
  print 'frq: ', oscbest.get_frq();
  print 'time up: ',oscbest.get_tim_up();
  print 'time down: ',oscbest.get_tim_down();
  print 'parameters', oscbest.get_rapams();
  print ""



  figure(1)
  clf()
  fnamet = time()
  t = oscbest.get_times();
  v1 =  oscbest.get_vals();
  v2 = WAVE.get_vals();
  t2 = range(len(WAVE.get_vals()));
  plot(t, v1, linewidth=1.0, color = 'blue')
  plot(t2, v2, linewidth=1.0, color = 'red')
  xlabel('time')
  ylabel('vals')
  xticks(range(0,TIME,int(TIME*0.1)))
  title(str(oscbest))
  grid(True)
  savefig('out/'+str(fnamet)+'rap'+str(STEP)+'.png',format='png')      

  figure(2)
  clf()
  plot(range(len(POP_SIZES)), POP_SIZES, linewidth=1.0, color = 'blue')
  savefig('out/'+str(fnamet)+'pop'+str(STEP)+'.png',format='png')

  if STEP == 0:
    if math.fabs(oscbest.get_amp()-WAVE.get_amp())/WAVE.get_amp() <= EPS_AMP:
      STEP = 1
    else:
			print "Repeating STEP 0..."
  elif STEP == 1:
    if math.fabs(oscbest.get_frq()-WAVE.get_frq())/WAVE.get_frq() <= EPS_FRQ:
      STEP = 2
      POP_MAX = POP_MAX*4
      POP_LIFETIME = POP_LIFETIME*2
    else:
			'''Ponovi step 1'''
			print "Repeating STEP 1..."
			#SAMPLING_EXP        = SAMPLING_EXP + 1
  elif STEP == 2:
    a = math.fabs(oscbest.get_amp()-WAVE.get_amp())/WAVE.get_amp() <= EPS_AMP
    b = math.fabs(oscbest.get_frq()-WAVE.get_frq())/WAVE.get_frq() <= EPS_FRQ
    c = math.fabs(oscbest.get_tim_up()-oscbest.get_tim_down())     <= EPS_TIM

    if a and b and c:
      STEP = 3
    else:
      print "Repeating STEP 2..."

  else:
      break;      

  rarray = deepcopy(oscbest.get_rapams())

print "End."
print rarray
